Spline plotting and SSANOVA
Stefano Coretta
Data import
Before importing the data, we need to specify the column names for the data set header, since the raw data does not have a header. The number of splines is always 42 (saved as num.splines). first.columns is a factor containing the names of columns that are not the splines coordinates columns. Finally, we can concatenate first.columns with the names of the splines coordinates which is generated by the paste0 function in the format X_1, Y_1, X_2, Y_2, X_n, ... (the underscore _ will be useful for separating the axis from the fan number).
num.splines <- 42
first.columns <- c(
"speaker",
"seconds",
"rec.date",
"prompt",
"label",
"TT.displacement.sm",
"TT.velocity",
"TT.velocity.abs",
"TD.displacement.sm",
"TD.velocity",
"TD.velocity.abs"
)
columns <- c(first.columns,
paste0(rep(c("X", "Y"), num.splines),
"_",
rep(1:num.splines, each = 2)
)
)We can now read in the file.
splines.raw <- list.files(path = "./pilot/results",
pattern = "*-aaa.csv",
full.names = TRUE) %>%
map_df(~read_tsv(., col_names = columns, na = "*"))## Parsed with column specification:
## cols(
## .default = col_double(),
## speaker = col_character(),
## rec.date = col_character(),
## prompt = col_character(),
## label = col_character()
## )## See spec(...) for full column specifications.## Parsed with column specification:
## cols(
## .default = col_double(),
## speaker = col_character(),
## rec.date = col_character(),
## prompt = col_character(),
## label = col_character()
## )## See spec(...) for full column specifications.## Parsed with column specification:
## cols(
## .default = col_double(),
## speaker = col_character(),
## rec.date = col_character(),
## prompt = col_character(),
## label = col_character()
## )## See spec(...) for full column specifications.## Parsed with column specification:
## cols(
## .default = col_double(),
## speaker = col_character(),
## rec.date = col_character(),
## prompt = col_character(),
## label = col_character()
## )## See spec(...) for full column specifications.rm(num.splines, first.columns, columns)
languages <- read_csv("./pilot/stimuli/languages.csv")## Parsed with column specification:
## cols(
## speaker = col_character(),
## language = col_character()
## )words <- read_csv("./pilot/stimuli/nonce.csv")## Parsed with column specification:
## cols(
## item = col_integer(),
## word = col_character(),
## ipa = col_character(),
## c1 = col_character(),
## c1phonation = col_character(),
## vowel = col_character(),
## anteropost = col_character(),
## height = col_character(),
## c2 = col_character(),
## c2phonation = col_character(),
## c2place = col_character(),
## language = col_character()
## )The following code applies tidy formatting to the data frame. It uses functions from the tidyr library.
splines <- splines.raw %>%
gather(spline, coordinate, matches("[XY]_")) %>%
separate(spline, c("axis", "fan"), convert = TRUE) %>%
spread(axis, coordinate) %>%
mutate(word = word(prompt, 2)) %>%
left_join(y = languages) %>%
left_join(y = words) %>%
mutate_if(is.character, as.factor)## Joining, by = "speaker"## Joining, by = c("word", "language")splines.it <- filter(splines, language == "italian")
splines.pl <- filter(splines, language == "polish")
splines.cor <- filter(splines, c2place == "coronal")
splines.vel <- filter(splines, c2place == "velar")Some plotting
We can finally plot splines.
filter(splines, label == "max_TD") %>%
ggplot(aes(x = X, y = Y, colour = c2phonation, group = word)) +
geom_smooth(method = "loess", se = FALSE) +
coord_fixed(ratio = 1) +
facet_grid(vowel ~ language)filter(splines, label == "max_TT") %>%
ggplot(aes(x = X, y = Y, colour = c2phonation, group = word)) +
geom_smooth(method = "loess", se = FALSE) +
coord_fixed(ratio = 1) +
facet_grid(vowel ~ language)SSANOVA
ssanovaTongue <- function(splines) {
model <- ssanova(Y ~ c2phonation + X + c2phonation:X, data = splines)
predicted <- expand.grid(X = seq(min(splines$X, na.rm = TRUE),
max(splines$X, na.rm = TRUE),
length = 100),
Y = seq(min(splines$Y, na.rm = TRUE),
max(splines$Y, na.rm = TRUE),
length = 100),
c2phonation = c("voiced", "voiceless")
)
predicted$splines.fit <- predict(model, predicted, se = T)$fit
predicted$splines.SE <- predict(model, predicted, se = T)$se.fit
return(predicted)
}splines.cor.a <- filter(splines.it, vowel == "a", label == "max_TT")
ssanova.cor.a <- ssanovaTongue(splines.cor.a)
splines.cor.o <- filter(splines.it, vowel == "o", label == "max_TT")
ssanova.cor.o <- ssanovaTongue(splines.cor.a)## Warning in chol.default(wk1, pivot = TRUE): the matrix is either rank-
## deficient or indefinitesplines.cor.u <- filter(splines.it, vowel == "u", label == "max_TT")
ssanova.cor.u <- ssanovaTongue(splines.cor.u)
splines.vel.a <- filter(splines.it, vowel == "a", label == "max_TD")
ssanova.vel.a <- ssanovaTongue(splines.vel.a)
splines.vel.o <- filter(splines.it, vowel == "o", label == "max_TD")
ssanova.vel.o <- ssanovaTongue(splines.vel.a)
splines.vel.u <- filter(splines.it, vowel == "u", label == "max_TD")
ssanova.vel.u <- ssanovaTongue(splines.vel.u)ggplot(ssanova.cor.a, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.cor.o, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.cor.u, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.a, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.o, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.u, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
splines.cor.a <- filter(splines.pl, vowel == "a", label == "max_TT")
ssanova.cor.a <- ssanovaTongue(splines.cor.a)
splines.cor.o <- filter(splines.pl, vowel == "o", label == "max_TT")
ssanova.cor.o <- ssanovaTongue(splines.cor.a)
splines.cor.u <- filter(splines.pl, vowel == "u", label == "max_TT")
ssanova.cor.u <- ssanovaTongue(splines.cor.u)## Warning in chol.default(wk1, pivot = TRUE): the matrix is either rank-
## deficient or indefinitesplines.vel.a <- filter(splines.pl, vowel == "a", label == "max_TD")
ssanova.vel.a <- ssanovaTongue(splines.vel.a)## Warning in chol.default(wk1, pivot = TRUE): the matrix is either rank-
## deficient or indefinite
## Warning in chol.default(wk1, pivot = TRUE): the matrix is either rank-
## deficient or indefinite
## Warning in chol.default(wk1, pivot = TRUE): the matrix is either rank-
## deficient or indefinitesplines.vel.o <- filter(splines.pl, vowel == "o", label == "max_TD")
ssanova.vel.o <- ssanovaTongue(splines.vel.a)
splines.vel.u <- filter(splines.pl, vowel == "u", label == "max_TD")
ssanova.vel.u <- ssanovaTongue(splines.vel.u)ggplot(ssanova.cor.a, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.cor.o, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.cor.u, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.a, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.o, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
ggplot(ssanova.vel.u, aes(x = X, colour = c2phonation)) +
geom_line(aes(y = splines.fit)) +
geom_ribbon(aes(ymin = splines.fit - (1.96 * splines.SE),
ymax = splines.fit + (1.96 * splines.SE),
fill = c2phonation
),
alpha = 0.5,
color = "NA"
)
Consonantal gestures: NONS, maximum, NOFS
Now we can create a separate data frame where the observasional unit is each word and the variables are the consonantal getures (target, maximum closure, release).
cons.gestures <- select(splines.raw, speaker:label) %>%
mutate(word = word(prompt, 2)) %>%
left_join(y = languages) %>%
left_join(y = words) %>%
separate(label, c("gesture", "tongue.area")) %>%
spread(gesture, seconds) %>%
mutate(nucleus = NOFF - NONS, deceleration = max - peak1) %>%
mutate_if(is.character, as.factor)## Joining, by = "speaker"## Joining, by = c("word", "language")gestures.it <- filter(cons.gestures, language == "italian")
gestures.pl <- filter(cons.gestures, language == "polish")Let’s plot the nucleus duration (NOFF - NONS) as a function of place of articulation, voicing of C2 (our target consonant in C1VC2V words), and language.
filter(cons.gestures, nucleus < 0.05) %>%
ggplot(aes(c2phonation, nucleus)) +
facet_grid(.~language + c2place) +
geom_violin() +
geom_boxplot(width=0.1) +
theme(strip.background = element_rect(fill = "black"),
strip.text = element_text(color = "white", face = "bold")
)
nucleus.model.it <- lmer(
nucleus ~
c2phonation *
vowel *
c2place +
(1|word) +
(1|speaker),
data = gestures.it
)
summary(nucleus.model.it)## Linear mixed model fit by REML ['lmerMod']
## Formula:
## nucleus ~ c2phonation * vowel * c2place + (1 | word) + (1 | speaker)
## Data: gestures.it
##
## REML criterion at convergence: -975.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.0662 -0.3538 -0.1067 0.1207 11.0059
##
## Random effects:
## Groups Name Variance Std.Dev.
## word (Intercept) 5.047e-06 0.002246
## speaker (Intercept) 8.853e-05 0.009409
## Residual 2.440e-04 0.015620
## Number of obs: 197, groups: word, 24; speaker, 2
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 0.019641 0.007609 2.581
## c2phonationvoiceless 0.008140 0.005269 1.545
## vowelo 0.001214 0.005323 0.228
## vowelu -0.002072 0.005523 -0.375
## c2placevelar -0.004562 0.006282 -0.726
## c2phonationvoiceless:vowelo -0.002168 0.007610 -0.285
## c2phonationvoiceless:vowelu -0.004056 0.007795 -0.520
## c2phonationvoiceless:c2placevelar -0.002087 0.008681 -0.240
## vowelo:c2placevelar 0.004645 0.008355 0.556
## vowelu:c2placevelar 0.004498 0.009085 0.495
## c2phonationvoiceless:vowelo:c2placevelar -0.003864 0.011852 -0.326
## c2phonationvoiceless:vowelu:c2placevelar -0.002311 0.012474 -0.185
##
## Correlation of Fixed Effects:
## (Intr) c2phnt vowelo vowelu c2plcv
## c2phntnvcls -0.339
## vowelo -0.336 0.485
## vowelu -0.324 0.468 0.463
## c2placevelr -0.287 0.409 0.407 0.389
## c2phonationvoiceless:vowelo 0.235 -0.692 -0.700 -0.324 -0.282
## c2phonationvoiceless:vowelu 0.229 -0.676 -0.328 -0.709 -0.275
## c2phntnvc:2 0.206 -0.607 -0.295 -0.284 -0.705
## vowelo:c2plcvlr 0.216 -0.308 -0.637 -0.293 -0.747
## vowelu:c2plcvlr 0.197 -0.284 -0.281 -0.607 -0.674
## c2phonationvoiceless:vowelo:c2plcvlr -0.151 0.444 0.449 0.207 0.520
## c2phonationvoiceless:vowelu:c2plcvlr -0.143 0.423 0.205 0.443 0.487
## c2phonationvoiceless:vowelo
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu 0.468
## c2phntnvc:2 0.420
## vowelo:c2plcvlr 0.444
## vowelu:c2plcvlr 0.197
## c2phonationvoiceless:vowelo:c2plcvlr -0.642
## c2phonationvoiceless:vowelu:c2plcvlr -0.293
## c2phonationvoiceless:vowelu c2ph:2
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2 0.410
## vowelo:c2plcvlr 0.207 0.530
## vowelu:c2plcvlr 0.430 0.486
## c2phonationvoiceless:vowelo:c2plcvlr -0.300 -0.733
## c2phonationvoiceless:vowelu:c2plcvlr -0.625 -0.696
## vowelo:c2plcvlr vowelu:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr 0.507
## c2phonationvoiceless:vowelo:c2plcvlr -0.701 -0.356
## c2phonationvoiceless:vowelu:c2plcvlr -0.367 -0.728
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phonationvoiceless:vowelu:c2plcvlr 0.509mixed(
nucleus ~
c2phonation *
vowel *
c2place +
(1|word) +
(1|speaker),
data = gestures.it
)## Fitting 8 (g)lmer() models:
## [........]
## Obtaining 7 p-values:
## [.......]## Effect df F.scaling F p.value
## 1 c2phonation 1, 7.68 1.00 2.39 .16
## 2 vowel 2, 8.97 1.00 0.17 .84
## 3 c2place 1, 15.42 1.00 0.52 .48
## 4 c2phonation:vowel 2, 9.03 1.00 0.14 .87
## 5 c2phonation:c2place 1, 14.13 1.00 0.06 .81
## 6 vowel:c2place 2, 13.46 1.00 0.18 .83
## 7 c2phonation:vowel:c2place 2, 13.04 1.00 0.05 .95plot(allEffects(nucleus.model.it))
nucleus.model.pl <- lmer(
nucleus ~
c2phonation *
vowel *
c2place +
(1|word),
data = gestures.pl
)
summary(nucleus.model.pl)## Linear mixed model fit by REML ['lmerMod']
## Formula: nucleus ~ c2phonation * vowel * c2place + (1 | word)
## Data: gestures.pl
##
## REML criterion at convergence: -487.2
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.3499 -0.6925 -0.1397 0.4125 4.1763
##
## Random effects:
## Groups Name Variance Std.Dev.
## word (Intercept) 0.0001650 0.01285
## Residual 0.0009129 0.03021
## Number of obs: 136, groups: word, 12
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 3.569e-02 1.534e-02 2.327
## c2phonationvoiceless 5.951e-03 2.157e-02 0.276
## vowelo 1.058e-03 2.182e-02 0.048
## vowelu -3.192e-03 2.157e-02 -0.148
## c2placevelar -1.547e-02 2.240e-02 -0.691
## c2phonationvoiceless:vowelo 1.009e-02 3.061e-02 0.330
## c2phonationvoiceless:vowelu 1.478e-03 3.043e-02 0.049
## c2phonationvoiceless:c2placevelar -4.777e-05 3.180e-02 -0.002
## vowelo:c2placevelar -2.947e-03 3.137e-02 -0.094
## vowelu:c2placevelar 3.914e-02 3.239e-02 1.208
## c2phonationvoiceless:vowelo:c2placevelar 4.969e-04 4.440e-02 0.011
## c2phonationvoiceless:vowelu:c2placevelar -2.210e-02 4.533e-02 -0.488
##
## Correlation of Fixed Effects:
## (Intr) c2phnt vowelo vowelu c2plcv
## c2phntnvcls -0.711
## vowelo -0.703 0.500
## vowelu -0.711 0.505 0.500
## c2placevelr -0.685 0.487 0.481 0.487
## c2phonationvoiceless:vowelo 0.501 -0.705 -0.713 -0.356 -0.343
## c2phonationvoiceless:vowelu 0.504 -0.709 -0.354 -0.709 -0.345
## c2phntnvc:2 0.482 -0.679 -0.339 -0.343 -0.704
## vowelo:c2plcvlr 0.489 -0.348 -0.696 -0.348 -0.714
## vowelu:c2plcvlr 0.474 -0.337 -0.333 -0.666 -0.692
## c2phonationvoiceless:vowelo:c2plcvlr -0.345 0.486 0.492 0.246 0.504
## c2phonationvoiceless:vowelu:c2plcvlr -0.338 0.476 0.238 0.476 0.494
## c2phonationvoiceless:vowelo
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu 0.500
## c2phntnvc:2 0.478
## vowelo:c2plcvlr 0.496
## vowelu:c2plcvlr 0.237
## c2phonationvoiceless:vowelo:c2plcvlr -0.689
## c2phonationvoiceless:vowelu:c2plcvlr -0.335
## c2phonationvoiceless:vowelu c2ph:2
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2 0.481
## vowelo:c2plcvlr 0.246 0.503
## vowelu:c2plcvlr 0.472 0.487
## c2phonationvoiceless:vowelo:c2plcvlr -0.345 -0.716
## c2phonationvoiceless:vowelu:c2plcvlr -0.671 -0.701
## vowelo:c2plcvlr vowelu:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr 0.494
## c2phonationvoiceless:vowelo:c2plcvlr -0.706 -0.349
## c2phonationvoiceless:vowelu:c2plcvlr -0.353 -0.714
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phonationvoiceless:vowelu:c2plcvlr 0.502mixed(
nucleus ~
c2phonation *
vowel *
c2place +
(1|word),
data = gestures.pl
)## Fitting 8 (g)lmer() models:
## [........]
## Obtaining 7 p-values:
## [.## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced## Warning in pf(Fstat, df1 = q, df2 = df2, lower.tail = FALSE): NaNs produced## ..## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced
## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced## ..## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced
## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced## .## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced
## Warning in pf(FstatU, df1 = q, df2 = df2, lower.tail = FALSE): NaNs
## produced## .]## Effect df F.scaling F p.value
## 1 c2phonation 1, 0 1.06 0.12 >.99
## 2 vowel 2, -0.00 -23053676236.74 <NA> <NA>
## 3 c2place 1, 0 0.98 -1.40 >.99
## 4 c2phonation:vowel 2, -0.00 -6029207709.50 <NA> <NA>
## 5 c2phonation:c2place 1, 0 0.97 0.00 >.99
## 6 vowel:c2place 2, -0.00 -158463611124.93 <NA> <NA>
## 7 c2phonation:vowel:c2place 2, -0.00 -66729713290.34 <NA> <NA>plot(allEffects(nucleus.model.pl))
ggplot(cons.gestures, aes(c2phonation, deceleration)) +
facet_grid(.~language + c2place) +
geom_violin() +
geom_boxplot(width=0.1) +
theme(strip.background = element_rect(fill = "black"),
strip.text = element_text(color = "white", face = "bold")
)## Warning: Removed 19 rows containing non-finite values (stat_ydensity).## Warning: Removed 19 rows containing non-finite values (stat_boxplot).
deceleration.model.it <- lmer(
deceleration ~
c2phonation *
vowel *
c2place +
(1|word) +
(1|speaker),
data = gestures.it
)
summary(deceleration.model.it)## Linear mixed model fit by REML ['lmerMod']
## Formula: deceleration ~ c2phonation * vowel * c2place + (1 | word) + (1 |
## speaker)
## Data: gestures.it
##
## REML criterion at convergence: -967.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.0341 -0.2432 0.1033 0.4218 2.9206
##
## Random effects:
## Groups Name Variance Std.Dev.
## word (Intercept) 8.083e-05 0.008991
## speaker (Intercept) 3.173e-04 0.017813
## Residual 8.202e-04 0.028639
## Number of obs: 250, groups: word, 24; speaker, 2
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 6.527e-02 1.537e-02 4.245
## c2phonationvoiceless 2.273e-05 1.247e-02 0.002
## vowelo 4.262e-03 1.254e-02 0.340
## vowelu -8.988e-03 1.272e-02 -0.707
## c2placevelar 1.155e-02 1.247e-02 0.926
## c2phonationvoiceless:vowelo 2.113e-03 1.773e-02 0.119
## c2phonationvoiceless:vowelu 9.166e-03 1.798e-02 0.510
## c2phonationvoiceless:c2placevelar 4.182e-03 1.768e-02 0.237
## vowelo:c2placevelar -1.909e-03 1.773e-02 -0.108
## vowelu:c2placevelar 1.972e-03 1.799e-02 0.110
## c2phonationvoiceless:vowelo:c2placevelar -5.694e-03 2.508e-02 -0.227
## c2phonationvoiceless:vowelu:c2placevelar -2.341e-02 2.538e-02 -0.922
##
## Correlation of Fixed Effects:
## (Intr) c2phnt vowelo vowelu c2plcv
## c2phntnvcls -0.405
## vowelo -0.403 0.497
## vowelu -0.397 0.490 0.487
## c2placevelr -0.405 0.500 0.497 0.490
## c2phonationvoiceless:vowelo 0.285 -0.703 -0.707 -0.344 -0.352
## c2phonationvoiceless:vowelu 0.281 -0.693 -0.345 -0.708 -0.347
## c2phntnvc:2 0.286 -0.705 -0.350 -0.345 -0.705
## vowelo:c2plcvlr 0.285 -0.351 -0.707 -0.345 -0.703
## vowelu:c2plcvlr 0.281 -0.346 -0.345 -0.707 -0.693
## c2phonationvoiceless:vowelo:c2plcvlr -0.201 0.497 0.500 0.244 0.497
## c2phonationvoiceless:vowelu:c2plcvlr -0.199 0.491 0.244 0.501 0.491
## c2phonationvoiceless:vowelo
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu 0.487
## c2phntnvc:2 0.496
## vowelo:c2plcvlr 0.500
## vowelu:c2plcvlr 0.244
## c2phonationvoiceless:vowelo:c2plcvlr -0.707
## c2phonationvoiceless:vowelu:c2plcvlr -0.345
## c2phonationvoiceless:vowelu c2ph:2
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2 0.489
## vowelo:c2plcvlr 0.244 0.496
## vowelu:c2plcvlr 0.500 0.489
## c2phonationvoiceless:vowelo:c2plcvlr -0.345 -0.705
## c2phonationvoiceless:vowelu:c2plcvlr -0.708 -0.697
## vowelo:c2plcvlr vowelu:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr 0.487
## c2phonationvoiceless:vowelo:c2plcvlr -0.707 -0.345
## c2phonationvoiceless:vowelu:c2plcvlr -0.345 -0.709
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phntnvcls
## vowelo
## vowelu
## c2placevelr
## c2phonationvoiceless:vowelo
## c2phonationvoiceless:vowelu
## c2phntnvc:2
## vowelo:c2plcvlr
## vowelu:c2plcvlr
## c2phonationvoiceless:vowelo:c2plcvlr
## c2phonationvoiceless:vowelu:c2plcvlr 0.491mixed(
deceleration ~
c2phonation *
vowel *
c2place +
(1|word) +
(1|speaker),
data = gestures.it
)## Fitting 8 (g)lmer() models:
## [........]
## Obtaining 7 p-values:
## [.......]## Effect df F.scaling F p.value
## 1 c2phonation 1, 11.31 1.00 0.00 >.99
## 2 vowel 2, 12.05 1.00 0.56 .59
## 3 c2place 1, 11.31 1.00 0.86 .37
## 4 c2phonation:vowel 2, 12.05 1.00 0.14 .87
## 5 c2phonation:c2place 1, 11.44 1.00 0.06 .82
## 6 vowel:c2place 2, 12.05 1.00 0.02 .98
## 7 c2phonation:vowel:c2place 2, 11.95 1.00 0.46 .64Kinematics performance
I want to check how good the gesture function in AAA is performing. The following gives the number of missing values per speaker.
cons.gestures %>%
select(speaker, max:peak1) %>%
group_by(speaker) %>%
summarise_each(funs(sum(is.na(.)))) %>%
select(speaker, peak1, NONS, max, NOFF)## # A tibble: 4 × 5
## speaker peak1 NONS max NOFF
## <fctr> <int> <int> <int> <int>
## 1 ab03 1 9 1 3
## 2 ps02 2 9 3 11
## 3 sc01 6 10 3 46
## 4 sdc02 1 6 6 4This chunk instead reports the total number of values per speaker.
cons.gestures %>%
select(speaker, max:peak1) %>%
group_by(speaker) %>%
summarise_each(funs(n())) %>%
select(speaker, peak1, NONS, max, NOFF)## # A tibble: 4 × 5
## speaker peak1 NONS max NOFF
## <fctr> <int> <int> <int> <int>
## 1 ab03 71 71 71 71
## 2 ps02 96 96 96 96
## 3 sc01 143 143 143 143
## 4 sdc02 120 120 120 120By comparing the two tables, it seems that the function is now performing quite well, but there could still be false detections.
GAMM
Coronal
splines.cor.max <- splines.cor %>%
filter(label == "max_TT") %>%
mutate(c2phonation.ord = as.ordered(c2phonation),
language.ord = as.ordered(language),
vowel.ord = as.ordered(vowel))
contrasts(splines.cor.max$c2phonation.ord) <- "contr.treatment"
contrasts(splines.cor.max$language.ord) <- "contr.treatment"
contrasts(splines.cor.max$vowel.ord) <- "contr.treatment"
coronal.gam <- bam(
Y ~
c2phonation.ord +
vowel.ord +
s(X, bs = "cr") +
s(X, by = c2phonation.ord, bs = "cr") +
s(X, by = vowel.ord, bs = "cr") +
s(rec.date, bs = "re") +
s(speaker, bs = "re"),
data = splines.cor.max,
method = "ML"
)
summary(coronal.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## Y ~ c2phonation.ord + vowel.ord + s(X, bs = "cr") + s(X, by = c2phonation.ord,
## bs = "cr") + s(X, by = vowel.ord, bs = "cr") + s(rec.date,
## bs = "re") + s(speaker, bs = "re")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -5.0345 0.3369 -14.942 < 2e-16 ***
## c2phonation.ordvoiceless 1.1511 0.2124 5.418 6.19e-08 ***
## vowel.ordo 2.0993 0.2585 8.121 5.29e-16 ***
## vowel.ordu 1.5753 0.2647 5.952 2.76e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(X) 8.751 8.956 787.234 < 2e-16 ***
## s(X):c2phonation.ordvoiceless 2.960 3.707 12.549 2.89e-09 ***
## s(X):vowel.ordo 4.601 5.627 12.267 9.73e-13 ***
## s(X):vowel.ordu 8.885 8.989 70.198 < 2e-16 ***
## s(rec.date) 84.337 201.000 0.734 9.00e-10 ***
## s(speaker) 2.528 3.000 17.326 1.48e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.745 Deviance explained = 74.8%
## -ML = 29627 Scale est. = 55.555 n = 8610coronal.gam.null <- bam(
Y ~
vowel.ord +
s(X, bs = "cr") +
s(rec.date, bs = "re") +
s(speaker, bs = "re"),
data = splines.cor.max,
method = "ML"
)
compareML(coronal.gam, coronal.gam.null)## coronal.gam: Y ~ c2phonation.ord + vowel.ord + s(X, bs = "cr") + s(X, by = c2phonation.ord,
## bs = "cr") + s(X, by = vowel.ord, bs = "cr") + s(rec.date,
## bs = "re") + s(speaker, bs = "re")
##
## coronal.gam.null: Y ~ vowel.ord + s(X, bs = "cr") + s(rec.date, bs = "re") + s(speaker,
## bs = "re")
##
## Chi-square test of ML scores
## -----
## Model Score Edf Chisq Df p.value Sig.
## 1 coronal.gam.null 30041.82 7
## 2 coronal.gam 29626.85 14 414.972 7.000 < 2e-16 ***
##
## AIC difference: -863.08, model coronal.gam has lower AIC.plot_smooth(coronal.gam, view="X", plot_all="c2phonation.ord", rug=F)## Summary:
## * c2phonation.ord : factor; set to the value(s): voiced, voiceless.
## * vowel.ord : factor; set to the value(s): a.
## * X : numeric predictor; with 30 values ranging from -53.195100 to 70.164600.
## * rec.date : factor; set to the value(s): 07/02/2017 16:21:26.
## * speaker : factor; set to the value(s): sc01.
plot_smooth(coronal.gam, view="X", plot_all="vowel.ord", rug=F)## Summary:
## * c2phonation.ord : factor; set to the value(s): voiceless.
## * vowel.ord : factor; set to the value(s): a, o, u.
## * X : numeric predictor; with 30 values ranging from -53.195100 to 70.164600.
## * rec.date : factor; set to the value(s): 07/02/2017 16:21:26.
## * speaker : factor; set to the value(s): sc01.
plot_diff(coronal.gam, view="X", comp=list(c2phonation.ord=c("voiceless","voiced")))## Summary:
## * vowel.ord : factor; set to the value(s): a.
## * X : numeric predictor; with 100 values ranging from -53.195100 to 70.164600.
## * rec.date : factor; set to the value(s): 07/02/2017 16:21:26.
## * speaker : factor; set to the value(s): sc01.
##
## X window(s) of significant difference(s):
## -53.195100 - -3.352797plot_diff(coronal.gam, view="X", comp=list(vowel.ord=c("o", "a")))## Summary:
## * c2phonation.ord : factor; set to the value(s): voiceless.
## * X : numeric predictor; with 100 values ranging from -53.195100 to 70.164600.
## * rec.date : factor; set to the value(s): 07/02/2017 16:21:26.
## * speaker : factor; set to the value(s): sc01.
##
## X window(s) of significant difference(s):
## -48.210870 - 4.123548
## 41.505276 - 60.196139acf_plot(resid(coronal.gam), split_by=list(splines.cor.max$rec.date))
Velar
splines.vel <- splines.vel %>%
mutate(c2phonation.ord = as.ordered(c2phonation))
contrasts(splines.vel$c2phonation.ord) <- "contr.treatment"
velar.gam <- bam(
Y ~
c2phonation.ord +
s(X, bs = "cr") +
s(X, by = c2phonation.ord, bs = "cr"),
data = splines.vel,
method = "ML"
)
summary(velar.gam)##
## Family: gaussian
## Link function: identity
##
## Formula:
## Y ~ c2phonation.ord + s(X, bs = "cr") + s(X, by = c2phonation.ord,
## bs = "cr")
##
## Parametric coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -4.17789 0.08214 -50.86 < 2e-16 ***
## c2phonation.ordvoiceless 0.72134 0.11559 6.24 4.42e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df F p-value
## s(X) 8.946 8.994 2321.441 < 2e-16 ***
## s(X):c2phonation.ordvoiceless 5.680 6.764 8.209 1.25e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.556 Deviance explained = 55.6%
## -ML = 1.2089e+05 Scale est. = 107.25 n = 32172velar.gam.null <- bam(
Y ~
c2phonation.ord +
s(X, bs = "cr"),
data = splines.vel,
method = "ML"
)
compareML(velar.gam, velar.gam.null)## velar.gam: Y ~ c2phonation.ord + s(X, bs = "cr") + s(X, by = c2phonation.ord,
## bs = "cr")
##
## velar.gam.null: Y ~ c2phonation.ord + s(X, bs = "cr")
##
## Chi-square test of ML scores
## -----
## Model Score Edf Chisq Df p.value Sig.
## 1 velar.gam.null 120911.0 4
## 2 velar.gam 120890.6 6 20.366 2.000 1.429e-09 ***
##
## AIC difference: -49.19, model velar.gam has lower AIC.plot_smooth(velar.gam, view="X", plot_all="c2phonation.ord", rug=F)## Summary:
## * c2phonation.ord : factor; set to the value(s): voiced, voiceless.
## * X : numeric predictor; with 30 values ranging from -56.350600 to 70.517200.
plot_diff(velar.gam, view="X", comp=list(c2phonation.ord=c("voiceless","voiced")))## Summary:
## * X : numeric predictor; with 100 values ranging from -56.350600 to 70.517200.
##
## X window(s) of significant difference(s):
## -43.535671 - -20.468798
## -10.216855 - 5.161061
## 16.694497 - 58.983764